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Abstract 

This paper investigates travelling wave solutions of the FitzHugh-Nagumo equation from the view- 
point of fast-slow dynamical systems. These solutions are homoclinic orbits of a three dimensional vector 
field depending upon system parameters of the FitzHugh-Nagumo model and the wave speed. Champ- 
neys et al. [A.R. Champneys, V. Kirk, E. Knobloch, B.E. Oldeman, and J. Sneyd, When Shilnikov meets 
Hopf in excitable systems, SIAM Journal of Applied Dynamical Systems, 6(4), 2007] observed sharp turns 
in the curves of homoclinic bifurcations in a two dimensional parameter space. This paper demonstrates 
numerically that these turns are located close to the intersection of two curves in the parameter space 
that locate non-transversal intersections of invariant manifolds of the three dimensional vector field. 
The relevant invariant manifolds in phase space are visualized. A geometrical model inspired by the 
numerical studies displays the sharp turns of the homoclinic bifurcations curves and yields quantitative 
predictions about multi-pulse and homoclinic orbits and periodic orbits that have not been resolved 
in the FitzHugh-Nagumo model. Further observations address the existence of canard explosions and 
mixed-mode oscillations. 

1 Introduction 

This paper investigates the three dimensional FitzHugh-Nagumo vector field defined by: 
ex\ = x 2 

ex 2 = ■jr(sx2-xt(xi-l)(a-xi)+y-p)=:—(sx 2 -f(x 1 )+y-p) (1) 

y = ~{xi-y) 
s 

where p, s, A, a and e are parameters. Our analysis views equations ([!]) as a fast-slow system with two 
fast variables and one slow variable. The dynamics of system ([I]) were studied extensively by Champneys et 
al. [1] with an emphasis on homoclinic orbits that represent travelling wave profiles of a partial differential 
equation pQ. Champneys et al. @] used numerical continuation implemented in AUTO [7] to analyze the 
bifurcations of ([!]) for e = 0.01 with varying p and s. As in their studies, we choose A = 5, a = 1/10 for 
the numerical investigations in this paper. The main structure of the bifurcation diagram is shown in Figure[T] 

Figure [T] shows a curve of Hopf bifurcations which is U-shaped and a curve of Shilnikov homoclinic bi- 
furcations which is C-shaped. Champneys et al. 14] observed that the C-curve is a closed curve which folds 
back onto itself before it reaches the U-curve, and they discussed bifurcations that can "terminate" a curve 
of homoclinic bifurcations. Their analysis does not take into account the multiple-time scales of the vector 
field 0. This paper demonstrates that fast-slow analysis of the homoclinic curve yields deeper insight into 
the events that occur at the sharp turn of the homoclinic curve. We shall focus on the turning point at the 
top end of the C-curvc and denote this region by /. 
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Figure 1: Bifurcation diagram for the FitzHugh-Nagumo equation 0. Shilnikov homoclinic bifurcations 
(solid red) and Hopf bifurcations (solid blue) are shown for e = 0.01. The dashed curves show the singular 
limit (e = 0) bifurcation curves for the homoclinic and Hopf bifurcations; see [15] and Section [5] for details 
on the singular limit part of the diagram. 



We regard e in the FitzHugh-Nagumo equation (JlJ as a small parameter. In [T5], we derived a singular 
bifurcation diagram which represents several important bifurcation curves in (p, s)-parameter space in the 
singular limit e = 0. The singular limits of the Hopf and homoclinic curves are shown in Figure [l] as dotted 
linesQ In the singular limit, there is no gap between the Hopf and homoclinic curves. We demonstrate below 
in Proposition 2.1 that a gap must appear for e > 0. The main point of this paper is that the termination 
point of the C-curve at the end of the gap is due to a fast-slow "bifurcation" where the two dimensional stable 
manifold of an equilibrium is tangent to the two dimensional unstable manifold of a one dimensional slow 
manifold^] Since the analysis of [4] does not explicitly consider slow manifolds of the system, this tangency 
does not appear in their list of possibilities for the termination of a C-curve. Note that the slow manifolds 
of the system are unique only up to "exponentially small" quantities of the form exp(— c/e),c > 0, so our 
analysis only identifies the termination point up to exponentially small values of the parameters. 

Fast-slow dynamical systems can be written in the form 

ex = ^=f(x,y,e) (2) 

where (x, y) G R m x R n and e is a small parameter < e < 1. The functions / :l'"xl"xl-) W n and 
g : R m x R" x K — > M. n are analytic in the systems studied in this paper. The variables x are fast and the 

1 In Section [5] we recall the precise meaning of the singular limit bifurcation from 115] and how they these bifurcations arise 
when e = 0. 

2 An analogous tangency plays an important role in the formation of mixed mode oscillations associated with singular Hopf 
bifurcations in fast-slow systems with one fast and two slow variables |12| . 
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variables y are slow. We can change (|2| from the slow time scale r to the fast time scale t = r/e, yielding 

x' = ~ = f(x,y,e) (3) 
y = = £9{x,y,e) 

In the singular limit e — > the system ([2]) becomes a differential-algebraic equation. The algebraic constraint 
defines the critical manifold: 

Co = {(x,y) eR m xl" : f(x,y,0) = 0} 

For a point p € Co we say that Co is normally hyperbolic at p if the all the eigenvalues of the m x m 
matrix D x f{p) have non-zero real parts. A normally hyperbolic subset of Co is an actual manifold and 
we can locally parametrize it by a function h(y) = x. This yields the slow subsystem (or reduced flow) 
y = g{h(y),y) defined on Co- Taking the singular limit e — > in ^ gives the fast subsystem (or layer 
equations) x' = f(x,y) with the slow variables y acting as parameters. Fenichel's Theorem [9] states that 
normally hyperbolic critical manifolds perturb to invariant slow manifolds C e . A slow manifold C e is 0(e) 
distance away from Co. The flow on the (locally) invariant manifold C e converges to the slow subsystem on 
the critical manifold as e — > 0. Slow manifolds are usually not unique for a fixed value of e = eo but lie at a 
distance O(e~ K / €0 ) away from each other for some K > 0; nevertheless we shall refer to "the slow manifold" 
for a fast-slow system with the possibility of an exponentially small error being understood. 

Section[2]discusses the fast-slow decomposition of the homoclinic orbits of the FitzHugh-Nagumo equation 
in the region /. This decomposition has been used to prove the existence of homoclinic orbits in the system 
for e sufficiently small [H [T71 HTJ US] , but previous work only applies to a situation where the equilibrium 
point for the homoclinic orbit is not close to a fold point. At a fold point the critical manifold of a fast-slow 
system is locally quadratic and not normally hyperbolic. This new aspect of the decomposition is key to 
understanding the sharp turn of the homoclinic curve. Section [3] presents a numerical study that highlights 
the geometric mechanism for the turning of the C-curve. We visualize relevant aspects of the phase portraits 
near the turns of the C-curve. In Section [4] we show that exponential contraction of the Shilnikov return 
map in the FitzHugh-Nagumo equation explains why n-homoclinic and n-periodic orbits are expected to 
be found at parameter values very close to a primary 1-homoclinic orbit. Section [5] presents two further 
observations. We identify where a canard explosion [26] occurs and we note the existence of two different 
types of mixed-mode oscillations in the system. 



2 Fast-Slow Decomposition of Homoclinic Orbits 

We introduce notation used in our earlier work |15j . The critical manifold of |T]) is given by: 

C = {(xi,x 2 ,y) e K 3 : x 2 = and y = f(x{) + p} 

It is normally hyperbolic away from the two fold points x\ t ± with x^- < a;i,+ which are found by solving 
f'(x\) = as the local minimum and maximum of the cubic /. Hence Co splits into three parts: 

C t = {xi < x lt _} n C , C m = {x-i_ < x x < xx,+} n C , C r = {x ly+ } n C 

We are mostly interested in the two branches C; and C r which are of saddle-type, i.e. points in C; and C r 
are saddle equilibria of the fast subsystem. The middle branch C m — {xi,±} consists of unstable foci for 
the fast subsystem. The slow manifolds provided by Fenichel's Theorem will be denoted by C; je and C r>e . 
The notation for the two-dimensional stable and unstable manifolds of C; j£ is W s (Ci ie ) and W u (C r ^) with 
similar notation for C r e ; the notation for the associated linear eigenspaces is e.g. E s (Ci ilL ). The full system 
([T]) has a unique equilibrium point which we denote by q. For (p, s) € / and e = 0.01 the dimensions of the 
stable and unstable manifolds are dim(W u (q)) — 1 and dim(W s (q)) = 2 with a complex conjugate pair of 
eigenvalues for the linearization at q. The equilibrium q is completely unstable inside the U-curve and the 
Hopf bifurcations we are interested in near / are all subcritical [HI |4] . 
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As e — > the Hopf bifurcation curve converges to a region in (p, s) parameter space bounded by two verti- 
cal lines p = p± and the segment {s = 0,p_ <p< P+}', see Figurc[I] The parameter values p± are precisely 
the values when the equilibrium point q coincides with the fold points Xi t ± |15j . This analysis gives one 
part of the singular limit bifurcation diagram showing what happens to the Hopf bifurcation curves for e = 0. 




Figure 2: Sketch of a homoclinic orbit to the unique equilibrium q. Fast (red) and slow (green) segments 
decompose the orbit into segments. 

When e is small, the homoclinic orbit in W u (q) n W s (q) can be partitioned into fast and slow segments. 
The singular limit of this fast-slow decomposition has four segments: a fast subsystem heteroclinic connection 
from g to C r , a slow segment on C r , a heteroclinic connection from C r to C; and a slow segment back to 
g on C;; see Figure [2] Existence proofs for the homoclinic orbits [5D1 [T71 [5] are based upon analysis of the 
transitions between these segments. Trajectories that remain close to a normally hyperbolic slow manifold 
must be "exponentially close" to the manifold except for short segments where the trajectory approaches the 
slow manifold along its stable manifold and departs along its unstable manifold. Existence of the homoclinic 
orbit depends upon how the four segments of its fast-slow decomposition fit together: 

(Fl) The one dimensional W u (q) approaches C r along its two dimensional stable manifold W s (C r ^). Inter- 
section of these manifolds cannot be transverse and occurs only for parameter values that lie along a 
curve in the (p, s) parameter plane. 

(51) The Exchange Lemma [TS] was developed to analyze the flow map for trajectories that approach C r ^ 
along its stable manifold and depart C r>e along its unstable manifold. 

(F2) The fast jump from a neighborhood of C Tt<L to a neighborhood of C;. e occurs along a transversal 
intersection of the two dimensional W s (Ci e ) and two dimensional W u (C re ). 

(52) The connection from C; j£ to q lies close to an intersection of the two dimensional W u (Ci j€ ) and the 
two dimensional W s {q). Previous analysis has dealt with parameter regions where the connection (S2) 
exists and is transversal, but it cannot persist up to the Hopf curve in the (p, s)-plane. 



Proposition 2.1. There exists a region in (p, s) -parameter space near the Hopf U-curve where no trajectories 
close to Ci >6 lie in W s (q). 

Proof. (Sketch) The Lyapunov coefficients of the Hopf bifurcations near / are positive [15], so the periodic 
orbits emanating from these bifurcations occur in the parameter region to the left of the Hopf curve. The 
periodic orbits are completely unstable. By calculating the eigenvalues of the linearization at the equilibrium 
we find that there is no Fold-Hopf bifurcation on the Hopf curve near /. Hence center manifold reduction 
implies that there will be a region of parameters near the Hopf curve where W s {q) is a topological disk 
whose boundary is the periodic orbit. Close enough to the Hopf curve, W s (q) and the periodic orbit lie at 
a finite distance from C;. £ and there is no connection from C;. £ to q. □ 
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This proposition implies that the parameter region in which there is a connection from Cj ;£ to q is bounded 
away from the Hopf curve. The next section shows that the boundary of this parameter region is very close 
to a curve along which there are tangential intersections of W u (Ci. e ) and W 8 {q). 

Remark: As e — >• 0, the C-curve converges to two lines (dashed red in Figure]!]) defined by homoclinic and 
heteroclinic orbits of the fast subsystem [15] . The horizontal segment of the C-curve to homoclinic orbits of 
the equilibrium point, and the sloped segment to heteroclinic orbits from the equilibrium point to the right 
branch of the critical manifold. Note that the C-curve terminates on the Hopf curve in the singular limit. 
The singular limit analysis does not explain the sharp turning of the C-curve for e > which is the focus of 
the next section. 

3 Interaction of Invariant Manifolds 

The slow manifold C'i_ t is normally hyperbolic away from the fold point X\— , with one attracting direction 
and one repelling direction. We recently introduced a method [H] for computing slow manifolds of saddle 
type. This algorithm is used here to help determine whether there are connecting orbits from a neighborhood 
of Cz i£ to the equilibrium point q. Our numerical strategy for finding connecting orbits has three steps: 

1. Choose the cross section 

S .09 = {{xi,x 2 ,y) eR 3 : y = 0.09} 

transverse to Ci ;C , 

2. Compute intersections of trajectories in W s (q) with S .og. These points are found either by backward 
integration from initial conditions that lie in a small disk D containing q in W s (q) or by solving a 
boundary value problem for trajectories that have one end in £0.09 and one end on the boundary of D. 

3. Compute the intersection pi € C; j£ f~l So. 09 with the algorithm described in Guckenheimer and Kuehn 
Q3| and determine the directions of the positive and negative eigenvectors of the Jacobian of the fast 
subsystem at p\. 



Figure |3] shows the result of these computations for e = 0.01, s — 1.37 and three values of p^\ The 
intersections of W s (q) with £0.09 lies close to W s (Ci^). 



e 


D=d(tangency,Hopf) 


1Q- 3 
10~ 4 


fa 1.07e 
ps l.OOe 
fa 0.98e 



Table 1: Euclidean distance in (p,s)-parameter space between the Hopf curve and the location of the tangency 
point between W s (q) and W u (C l>e ). 

Backward trajectories flowing along C;. £ converge to its stable manifold at a fast exponential rate. This 
fact also explains the observation that W s (q) H So. 09 makes a sharp turn. In Figure ^ a), it is apparent that 
the turn lies to the left of W u (Q, e ) H £ .09 and that W s {q) H W u {Ci t€ ) is non-empty. In Figure ^c) , the 
turn lies to the right of W u (Ci je ) D S .o9- We have also computed the distance from the Hopf curve of the 
parameters at which W s (q) and W u {Ci^) appear to have a tangential intersection for several different values 
of e; see Table IT] from which we observe that the distance is O(e). 



3 The second step above was carried out with two different initial value solvers, odel5s in Matlab I31| and dop853 |16| . and 
with AUTO [7] used as a boundary value solver. All three methods produced similar results. 
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(c) p = 0.06, s = 1.37 (d) Parameter space, region /. 



Figure 3: Figures (a)-(c) show the movement of the stable manifold W s (q) (cyan) with respect to E u (Ci^) 
(red) and E s (Ci. c ) (green) in phase space on the section y = 0.09 for e = 0.01. The parameter space 
diagram (d) shows the homoclinic C-curve (solid red), an extension of the C-curve of parameters where 
W u (q) Pi W s (C r>e ) is nonempty, a curve that marks the tangency of W s {q) to E u (Ci. e ) (blue) and a curve 
that marks a distance between C/ i£ and W s (q) (dashed blue) of 0.01 where the arrows indicate the direction 
in which the distance is bigger than 0.01. The solid black squares in (d) show the parameter values for 
(a)-(c). 
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In Figure |3[d) the C-curve of homoclinic bifurcations (solid red) has been computed using continuation 
in AUTO [7] as carried out by Champneys et al. [4, . Despite the fact that no homoclinic orbit exists in part 
of the region / it is possible to check whether the unstable manifold W u (q) reaches a small neighborhood of 
W s (C r}£ ). This idea has been used in a splitting algorithm |TS] to calculate where homoclinic orbits would 
occur if W s (q) would not move away from C; j£ as shown in Figures [3^a)-[3](d) . This yields the dashed red 
curve in Figure [3][d). On this curve we verified that W s {Ci^) and W u {C r<f ) still intersect transversally by 
computing those manifolds; see [IS] H3] for details. 

The blue curves in Figure |3]jd) have been obtained by measuring the Euclidean distances between W s (q) 
and Ci <e in the section So.09- Along the dashed blue curve the distance between C/ ie and W s (q) is 0.01. 
The arrows indicate the direction in which this distance increases. The solid blue curve marks a tangency 
of W s (q) with E u (Ci. t ). These calculations demonstrate that the sharp turn in the C-curve of homoclinic 
bifurcations occurs very close to the curve where there is a tangential intersection of W s (q) and W"(Cz >e ). 
Therefore, we state the following conjecture. 

Conjecture 3.1. The C-curve of homoclinic bifurcations of the FitzHugh-Nagumo system turns exponentially 
close to the boundary of the region where W u (Ci e ) f~l W s (q) is nonempty. 



Note that trajectory segments of types (Fl), (SI) and (F2) are still present along the dashed red curve 
in Figure |3jd). Only the last slow connection (S2) no longer exists. Existence proofs for homoclinic orbits 
that use Fenichel's Theorem for C; to conclude that trajectories entering a small neighborhood of C;. e must 
intersect W s (q) break down in this region. The equilibrium q has already moved past the fold point Xi t — in 
/ as seen from the singular bifurcation diagram in Figure [I] where the blue dashed vertical lines mark the 
parameter values where q passes through X\ ±. Therefore Fenichel's Theorem does not provide the required 
perturbation of Q e . Previous proofs [501 [HI 12] assumed p — and the connecting orbits of type (S2) do 
exist in this case. 

Shilnikov proved that there are chaotic invariant sets in the neighborhood of homoclinic orbits to a 
saddle-focus in three dimensional vector fields when the magnitude of the real eigenvalue is larger than the 
magnitude of the real part of the complex pair of eigenvalues |32j . The homoclinic orbits of the FitzHugh- 
Nagumo vector field satisfy this condition in the parameter region /. Therefore, we expect to find many 
periodic orbits close to the homoclinic orbits and parameters in I with "multi-pulse" homoclinic orbits that 
have several jumps connecting the left and right branches of the slow manifold [8]. Without making use 
of concepts from fast-slow systems, Champneys et al. [4] described interactions of homoclinic and periodic 
orbits that can serve to terminate curves of homoclinic bifurcations. This provides an alternate perspective 
on identifying phenomena that occur near the sharp turn of the C-curve in I. AUTO can be used to locate 
families of periodic orbits that come close to a homoclinic orbit as their periods grow. 

Figure H shows several significant objects in phase space for parameters lying on the C-curve. The homo- 
clinic orbit and the two periodic orbits were calculated using AUTO. The periodic orbits were continued in 
p starting from a Hopf bifurcation for fixed s ~ 1.3254. Note that the periodic orbit undergoes several fold 
bifurcations 0], We show two of the periodic orbits arising at p = 0.05; see [J|. The trajectories in TU s (Cz je ) 
have been calculated using a mesh on C/ j£ and using backward integration at each mesh point and initial 
conditions in the linear approximation E s {Ci^). 

We observe from Figure 4] that part of W s {q) lies near as expected for (S2) to be satisfied. This is 
in contrast to the situation beyond the turning of the C-curve shown in Figure [5]for p= 0.06 and s = 1.38. 
We observe that W s (q) is bounded. Figure [5|a) shows two periodic orbits Pi and P2 obtained from a Hopf 
bifurcation continuation starting for s = 1.38 fixed. P2 is of large amplitude and is obtained after the first 
fold bifurcation occurred. Pi is of small amplitude and is completely unstable. A zoom near P\ in Figure 
[5jb) and a time series comparison of a trajectory in W s (q) and P\ in Figure [5jc) show that 

lim{p : p g W s (q) and p ^ q} = P x (4) 
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Figure 4: Phase space along the C-curve near its sharp turn: the parameter values e = 0.01, p = 0.05 and 
s m 1.3254 lie on the C-curve. The homoclinic orbit (red), two periodic orbits born in the subcritical Hopf 
(blue), Co (thin black), Ci j£ and C r<e (thick black) are shown. The manifold W s (q) (cyan) has been truncated 
at a fixed coordinate of y. Furthermore W s (C^ e ) (green) is separated by C; >e into two components shown 
here by dark green trajectories interacting with C m ^ and by light green trajectories that flow left from C\^. 

where lim Q U denotes the a-limit set of some set U C R m x W\ From ^ we can also conclude that there 
is no hcteroclinic connection from q to Pi and only a connection from Pi to q in a large part of the region I 
beyond the turning of the C-curve. Since Pi is completely unstable, there can be no heteroclinic connections 
from q to Pi. Therefore, double heteroclinic connections between a periodic orbit and q are restricted to 
periodic orbits that lie closer to the homoclinic orbit than P±. These can be expected to exist for parameter 
values near the end of the C-curve in accord with the conjecture of Champneys et al. 4| and the "Shilnikov"- 
modcl presented in the next section. 

Remark: The recent manuscript [3] extends the results of [1] that motivated this paper. A partial 
unfolding of a heteroclinic cycle between a hyperbolic equilibrium point and a hyperbolic periodic orbit is 
developed in [3] . Champneys et al. call this codimension two bifurcation an EP lt-cycle and the point 
where it occurs in a two dimensional parameter space an EPlt-point. The manuscript [3] does not conclude 
whether the EP lt-scenario occurs in the FitzHugh-Nagumo equation. The relationship between the results 
of this paper and those of [3] have not yet been clarified. 

4 Homoclinic Bifurcations in Fast-Slow Systems 

It is evident from Figure [3] that the homoclinic orbits in the FitzHugh-Nagumo equation exist in a very thin 
region in (p, s)-parameter space along the C-curve. We develop a geometric model for homoclinic orbits that 
resemble those in the FitzHugh-Nagumo equation containing segments of types (SI), (Fl), (S2) and (F2). 
The model will be seen to be an exponentially distorted version of the Shilnikov model for a homoclinic orbit 
to a saddle- focus [13J . Throughout this section we assume that the parameters lie in a region / the region 
of the (p, s)-plane close to the upper turn of the C-curve. 

The return map of the Shilnikov model is constructed from two components: the flow map past an 
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Figure 5: The parameter values are e = 0.01, p — 0.06 and s = 1.38. For (a) we display two periodic orbits 
(blue), one with a single large excursion P 2 and one consisting of a small loop Pi. We also show q (red dot), 
trajectories in W s (Ci ie ) (green) and W s (q) (cyan). In (b) a zoom near q is shown and we made plotted a 
single trajectory 7 € W s (q) (cyan). The plot (c) gives a time series of this trajectory 7 in comparison to 
the periodic orbit P\. Note that the trajectories are computed backward in time, so the final points of the 
trajectories are on the left of the figure. A phase shift of time along the periodic orbit would bring the two 
time series closer. 
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equilibrium point, approximated by the flow map of a linear vector field, composed with a regular map that 
gives a "global return" of the unstable manifold of the equilibrium to its stable manifold [T3]. Place two 
cross-sections Si and S2 moderately close to the equilibrium point and model the flow map from Si to £2 
via the linearization of the vector field at the equilibrium. 




Figure 6: Sketch of the geometric model for the homoclinic bifurcations. Only parts of the sections for 
i = l,2 are shown. 



The degree one coefficient of the characteristic polynomial at the equilibrium has order O(e), so the imag- 
inary eigenvalues at the Hopf bifurcation point have magnitude (^(e 1 / 2 ). The real part of these eigenvalues 
scales linearly with the distance from the Hopf curve. Furthermore we note that the real eigenvalue of the 
equilibrium point remains bounded away from as e — > 0. 

Let ip(xi, X2, y) = (u,v,w) be a coordinate change near q so that tp(q) = and the vector field is in 
Jordan normal form up to higher order terms. We denote the sections obtained from the coordinate change 
into Jordan form coordinates by Si = V'(Si) and S2 = ^(£2); see FigureJHJ Then the vector field is 

v! — —j3u — av 

v = au — j3v + h.o.t. (5) 
w' = 7 

with a, (3, 7 positive. We can choose so that the cross-sections are S x = {u = 0, w > 0} and S 2 = {w = 1}. 
The flow map F12 : Si — > S2 of the (linear) vector field ^ without higher-order terms is given by 

Fi2(v,w) = uw^ 7 ^cos ( — ln(iu)J , sin ( — ln(u>)J J (6) 

Here j3 and a tend to as e — > 0. The domain for F12 is restricted to the interval v g [exp(— 2Tr/3/a), 1] 
bounded by two successive intersections of a trajectory in W s (0) with the cross-section u = 0. 

The global return map R : S2 —¥ Si of the FitzHugh-Nagumo system is obtained by following trajectories 
that have successive segments that are near W s (C r ^) (fast), CV j£ (slow), W v (C r<e ) fl W s (Ci^) (fast), C/ i£ 
(slow) and W u (Ci^) (fast). The Exchange Lemma [12] implies that the size of the domain of R in S2 is a 
strip whose width is exponentially small. As the parameter p is varied, we found numerically that the image 
of R has a point of quadratic tangency with W s (q) at a particular value of p. We also noted that W u (q) 
crosses W s (C rte ) as the parameter s varies [15]. Thus, we choose to model R by the map 

(w, v) = F 21 (u, v) = {av + A 2 - p 2 (u - Xif , p{u - Aj) + A 3 ) (7) 

for F21 where Ai represents the distance of W u (q) H S2 from the domain of F21, X2 represents how far the 
image of F 2 i extends in the direction normal to W s (q), A3 is the v coordinate of i*2i(Ai,Q) and p^ 1 , a are 
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0{e- K / e ) for suitable K > 0. Wc assume further that the domain of F21 is [Ai,Ai + p x [—1, 1]. Figure 
[7] depicts F21. With these choices, we observe two properties of the C-curve of homoclinic orbits in the 
geometric model: 

1. If cry + X 2 — p 2 {u — Ai) 2 is negative on the domain of F 2 i, then the image of F 21 is disjoint from the 
domain of F\ 2 and there are no recurrent orbits passing near the saddle point. Thus, recurrence implies 
that A2 > —a. 

2. If A2 > 0, then there are two values of Ai for which the saddle-point has a single pulse homoclinic orbit. 
These points occur for values of Ai for which the w-component of F 2 i(0, 0) vanishes: Ai = 1 1 A2 1 1//2 - 
The magnitude of these values of A2 is exponentially small. 



V > 










u 

*■ 






Ai 


W u {q) 







F< 



21 




(A 2 ,A 3 ) 



Figure 7: Sketch of the map F 21 : S2 — > Si- The coordinates are centered at W u (q) and the domain 

of F 2 i is in the thin rectangle at distance Ai from the origin. The image of this rectangle is the parabolic 
strip in Ei. 



When a vector field has a single pulse homoclinic orbit to a saddle-focus whose real eigenvalue has larger 
magnitude than the real part of the complex eigenvalues, Shilnikov |32j proved that a neighborhood of this 
homoclinic orbit contains chaotic invariant sets. This conclusion applies to our geometric model when it 
has a single pulse homoclinic orbit. Consequently, there will be a plethora of bifurcations that occur in the 
parameter interval A2 G [0,ct], creating the invariant sets as A2 decreases from a to 0. 

The numerical results in the previous section suggest that in the FitzHugh-Nagumo system, some of the 
periodic orbits in the invariant sets near the homoclinic orbit can be continued to the Hopf bifurcation of the 
equilibrium point. Note that saddle-node bifurcations that create periodic orbits in the invariant sets of the 
geometric model lie exponentially close to the curve A2 = that models tangency of W s (q) and W u (Ci te ) 
in the FitzHugh-Nagumo model. This observation explains why the right most curve of saddle-node bifur- 
cations in Figure 7 of Champneys et al. [4] lies close to the sharp turn of the C-curve. 

There will also be curves of heteroclinic orbits between the equilibrium point and periodic orbits close to 
the C-curve. At least some of these form codimension two EPlt bifurcations near the turn of the C-curve 
as discussed by Champneys et al. [I]. Thus, the tangency between W s (q) and W u (Ci^) implies that there 
are several types of bifurcation curves that pass exponentially close to the sharp turn of the C-curve in the 
FitzHugh-Nagumo model. Numerically, any of these can be used to approximately locate the sharp turn of 
the C-curve. 



5 Canards and Mixed Mode Oscillations 

This section reports two additional observations about the FitzHugh-Nagumo model resulting from our 
numerical investigations and analysis of the turning of the C-curve. 



11 



5.1 Canard Explosion 



The previous sections draw attention to the intersections of W s (q) and W u (Cj i£ ) as a necessary component 
for the existence of homoclinic orbits in the FitzHugh-Nagumo system. Canards for the backward flow of 
this system occur along intersections of W u (C^ e ) and C m ^. These intersections form where trajectories that 
track C; i£ have continuations that lie along C m ^ which has two unstable fast directions. We observed from 
Figures [3] and [5] that a completely unstable periodic orbit born in the Hopf bifurcation on the U-curve under- 
goes a canard explosion, increasing its amplitude to the size of a relaxation oscillation orbit upon decreasing 
p. This canard explosion happens very close to the intersections of W u (Ci :e ) and C m , e . 




0.0585 0.059 
P 




(a) (p : s)-space: Black circles correspond to two portraits in (b). 



(b) (x\ , y)-projection. 



Figure 8: The dashed green curve indicates where canard orbits start to occur along C mj£ . For values of 
p to the left of the dashed green curve we observe that orbits near the middle branch escape in backward 
time (upper panel in (b)). For values of p to the right of the dotted green curve trajectories near C m ^ stay 
bounded in backward time. 



To understand where this transition starts and ends we computed the middle branch C m>e of the slow 
manifold by integrating backwards from points between the fold points x\- and £1.+ starting close to C m ,o 
and determined which side of W M (Ci^) these trajectories came from. The results are shown in Figure [8] 
The dashed green curve divides the (p, s) plane into regions where the trajectory that flows into C m;€ lies to 
the left of W"(C| i£ ) and is unbounded from the region where the trajectory that flows into C m ^ lies to the 
right of W u (Ci ie ) and comes from the periodic orbit or another bounded invariant set. This boundary was 
found by computing trajectories starting on C' m .o backward in time. In backward time the middle branch of 
the slow manifold is attracting, so the trajectory first approaches C m ^ and then continues beyond its end 
when x\ decreases below X\ -. Figure |H] illustrates the difference in the behavior of these trajectories on the 
two sides of the dashed green curve. It shows that the parameters with canard orbits for the backward flow 
have smaller values of p than those for which W s (q) and W u (Ci^) have a tangential intersection. The turns 
of the C-curve do not occur at parameters where the backward flow has canards. 



5.2 Mixed-Mode Oscillations 

Mixed-mode oscillations (MMOs) have been observed in many fast-slow systems; see e.g. [28 | \29 \ \3U \ 112 ) . 
MMOs are periodic orbits which consist of sequences of small and large amplitude oscillations. The notation 
L s is used to indicate an MMO with L large and s small oscillations. 

The FitzHugh-Nagumo equation exhibits MMOs: the periodic orbits close to the homoclinic orbit 
make small oscillations near the equilibrium point in addition to large amplitude relaxation oscillations. A 
l 1 MMO is shown in Figure [9^ a)- (b). It was obtained by switching from the homoclinic C-curve to a nearby 



12 



O.Ur 
0.12 

0.1 - 
0.08 - 
0.06 - 
0.04 - 
0.02 - 



- 
-0.4 





(a) (a;i, j/)-space, p = 0.00266172. 



(b) Time series for (a). 



0.2 
0.18 
0.16 

o.u 

yO.12 
0.1 
0.08 
0.06 




0.2 0.4 

Xl 




(c) (xi,i/)-space, p = 0.0628718. 



(d) Time Series for (c). 



Figure 9: Some examples of mixed-mode oscillations in the FitzHugh-Nagumo equation. Fixed parameter 
values are e = 0.01 and s = 1. Note that the period of the orbits has been rescaled to 1 in (&) and (d). 
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curve of periodic orbits in a continuation framework. Note that the existence of multi-pulse homoclinic orbits 
near a Shilnikov homoclinic orbit [TO] H] implies that much more complicated patterns of MMOs also exist 
near the homoclinic C-curve. L s MMOs with very large L and s near the homoclinic C-curve are theoreti- 
cally possible although observing them will be very difficult due to the exponential contraction described in 
Section HI 

In addition to the MMOs induced by the Shilnikov bifurcation we also find MMOs which exist due to 
orbits containing canard segments near the completely unstable slow manifold G m ^. An example of a 4 1 
MMO is shown in Figure |9jc)-(d) obtained by continuation. In this case the small oscillations arise due 
to small excursions reminiscent to MMOs in three-time scale systems jT5J [53]. MMOs of type L 1 with 
L = 1, 2, 3, . . . , O(10 2 ) can easily be observed from continuation and we expect that L 1 MMOs exist for any 
L £ N. It is likely that these MMOs can be analyzed using a version of the FitzHugh-Nagumo equation 
containing 0(1), O(e) and 0(e 2 ) terms similar to the one introduced in [15] but we leave this analysis for 
future work. 

Figure[9]was obtained by varying p for fixed values of e — 0.01 and s — 1. Thus, varying a single parameter 
suffices to switch between MMOs whose small amplitude oscillations have a different character. In the first 
case, the small amplitude oscillations occur when the orbit comes close to a saddle focus rotating around 
its stable manifold, while in the second case, the trajectory never approaches the equilibrium and its small 
amplitude oscillations occur when the trajectory flows along the completely unstable slow manifold C m ^ e . 
Different types of MMOs seem to occur very frequently in single- and multi-parameter bifurcation problems; 
see [5] for a recent example. This contrasts with most work on the analysis of MMOs [551 US] that focuses 
on identifying the mechanism for generating MMOs in an example. The MMOs in the FitzHugh-Nagumo 
equation show that a fast-slow system with three or more variables can exhibit MMOs of different types and 
that one should not expect a priori that a single mechanism suffices to explain all the MMO dynamics. 
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